/*! \file functions.hpp
 * \brief routines which define functions
 */

#ifndef FUNCTIONS_HPP_
#define FUNCTIONS_HPP_

#include <math.h>
#include "evector.hpp"
#include "gauss.hpp"
#include "superelement.hpp"
#define PI (3.14159265358979323846)
const double sigma0 = 0;    //vd1
const double sigma1 = 1;

/*! \brief returns a value of Type T for each position vector
 * \param[in] x position vector of the evaluation point
 * \return value of type T assigned to x
 *
 * Template argument is the type of the computed values of the function (e.g. double)
 */
template<typename T>  T ExactSolution(double x)
{
  return sin(PI*x);
}

template<typename T>  T Function_p(double x)
{
  return 1.0;
  return x*x+1.0;
}

template<typename T>  T Function_q(double x)
{
  return 0.0;
  return x*x+1.0;
  return exp(x);
}

template<typename T>  T Function_f(double x)
{
	return 10;
	return x+2.0;
//	return x*(1.0-x)+2.0;
//	return (2.0+x*x)*sin(x);
  return (PI*PI+1.0)*sin(PI*x);
//  return sin(PI*x)-2.0*PI*x*cos(PI*x)+PI*PI*(1.0+x*x)*sin(PI*x);
  return exp(x)*sin(PI*x)-2.0*PI*x*cos(PI*x)+PI*PI*(1.0+x*x)*sin(PI*x);
}

template<typename T>  T Function_Bkk(double a, double b, double x)
{
	double h = b-a;
	double y = (b-a)*x + a;
	return h * (Function_p<double>(y)/pow(h, 2.0) + Function_q<double>(y)*pow(x, 2.0));
}
template<typename T>  T Function_Bk1k1(double a, double b, double x)
{
	double h = b-a;
	double y = (b-a)*x + a;

	return h * (Function_p<double>(y)/pow(h, 2.0) + Function_q<double>(y)*pow(1-x, 2.0));
}
template<typename T>  T Function_Bkk1(double a, double b, double x)
{
	double h = b-a;
	double y = (b-a)*x + a;

	return h * (-Function_p<double>(y)/pow(h, 2.0) + Function_q<double>(y)*x*(1-x));
}
template<typename T>  T Function_Fk(double a, double b, double x)
{
	double h = b-a;
	double y = (b-a)*x + a;
	return h * Function_f<double>(y)*x;
}

template<typename T>  T Function_Fk1(double a, double b, double x)
{
	double h = b-a;
	double y = (b-a)*x + a;
	return h * Function_f<double>(y)*(1-x);
}


//********************* define a basic function ************************
/*! \brief constant basis function defined on the interval [0,1]
 * \param[in] x position in interval [0,1]
 * \return value 1 
 */
template<typename T>  T shapeFunction0 (double x) {return x;}

/*! \brief linear shape function (1-x) defined on the interval [0,1]
 * \param[in] x position in interval [0,1]
 * \return value (1-x)
 */
template<typename T>  T shapeFunction1_0 (double x) {return (1.0-x);}

/*! \brief linear shape function (x) defined on the interval [0,1]
 * \param[in] x position in interval [0,1]
 * \return value (x)
 */

#endif // FUNCTIONS_HPP_
